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An improved method for the calculation of asymptotic phases in quantum 
mechanical scattering problems is developed. The method is based on the nu- 
merical integration of the radial wave equation by using a new modification 
of Gauss’s central difference formula. Applications are made to the computa- 
tion of some S-phases connected with the Yukawa potential, and a comparison 
with the results obtained by HutrHen according to his variational principle 
shows that his simplified modification of this principle gives phase values which 
are extremely accurate even when only a few variational parameters are used. 


I. Introduction 


The asymptotic phases 7 are of essential importance for the calculation of 
cross-sections in the quantum mechanical theory of scattering!, and in the 
literature various methods have been developed for their computation by means 
of numerical integration.’ * Since many of these methods are rather laborious 
to handle practically, we will here give a brief description of a method of 
numerical integration which is extremely simple and accurate. 


Il. A modification of Gauss’s formulas for numerical integration 


The radial wave equation of the scattering problem of a central field is an 
ordinary differential equation of the second order of the type 


ag — H(z, y) (1) 


with the initial condition y = yy, y’ = yo for «=a. We will use the theory 
of finite differences for the numerical integration of this equation. Let h be the 


* Now on leave at Department of Physics, Duke University, Durham, North Carolina, 

1H. Faxtn—J. Horrsmark, Zeits. f. Physik 45, 307 (1927). 

2 V. Rosgansky—W. WEeETzEL, Phys. Rev. 38, 1979 (1931); P. M. Morsz, Rev. Mod. Phys. 
4, 577 (1932); D. R. Harrren, R..pE L. Kronic, and H. Petersen, Physica 1, 895 (1934); 
R. G. Sacus and M. Gorprert-Mayer, Phys. Rev. 53, 991 (1938); and others. 
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interval, and let E be the step operator, defined by H/(x)=/(a +h). Then 


we will introduce Suepparp’s central difference operator 6 and mean-value — 


operator « by 


j= —F + pas. PE (2) 


sie 


with uw? = 1+ 62/4, and w(1 + 62/4)-?=1. Following SHepparp® we will use 
operator calculus and symbolic expansions. Taylor’s theorem gives immediately 
E=exp (hD), where D is the differentiation operator D=d/dz. Further we 
get from eq. (2) that hD—= 2 sinh! 6/2. The differential equation (1) can there- 
fore be solved by the inverse formula 


hg= {2 amh—! 6/2)-2 7 = 


(3) 
l 1 31 289 
eae (ee 2 | 4 6 fo se 
maar i? ie cs 60 480° 3628.800° 3 
and the derivative y’ can, if required, be obtained from 
h-ty’ = (1 + 62/4)-? {2 sinh-16/2}-1 F = a 
vs ll 191 2497 
pee, —1 ct ta ash 3 5 7 ~sh 
See gh! ont angi, (eee | 


Here the symbols 6-!F and 6-?F denote the first and second sums of F, re- 
spectively, ie. 6(6-1F) = F etc. The “integration constants” included in these 
sums will here be determined by the condition that the formulas under con- 
sideration should be valid even in the starting point, x=2). The equations 
(3) and (4) were both given by Gauss, and they have been of fundamental 
importance for the theory of numerical integration using finite differences. 
They are here derived only symbolically, but they can also be justified and 
their remainder-terms can be calculated by using Newton’s interpolation for- 
mula for unequal intervals with Cauchy’s remainder.® 

Gauss’s formula (3) can directly be employed for the calculation of asymp- 
totic phases®, but here we will instead use a slightly modified form developed 
by one of us (P.O.L.). The radial wave equation is linear, and we will therefore 
confine ourselves to considering only the simplified form 


dq? 
ree (x) y = f (2). (5) 


3 W. F. SHePpaRD, Proc. London Math. Soc. 31, 449 (1899). See also J. G. L. Micuert, 
J. Inst. Actuar. 72, 470 (1946); W. G. BickitEy, J. Math. Physics 27, 183 (1948). 

4 Encke, Berliner Astr. Jahrbuch 1837, ibid. 1867; CownLt1i, Greenwich Obs. 1909; C. SvéR- 
MER, Videnskapsselskapets skr. Oslo, Mat. natv. kl. 1913, no. 4, 10, and 14; J. Jackson, 
Mon. Not. Roy. Astr. Soc. 84, 602 (1924); D. R. Harrrer, Mem. Proc. Lit. Phil. Soc. Man- 
chester 77, 91 (1933). 

° H. P. Nietsen, Arkiv f. mat. fysik o. astr. 4, no. 21 (1908); E. J. Nystrém, Acta 
Soc. Fennicae 50, no, 13 (1926); J. F. Srerrensen, Skan. Aktuarietidskr. 7, 1 (1924); and 
others. 

® See, for instance, R. A. BuckineHAmM in Morr and Massry, Theory of Atomic Collisions, 
2nd ed. (Oxford 1949), p. 129. 
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This special case has also been treated by other authors.? Transforming the 
F-term in eq. (3) to the left member’, we obtain 


vi ts g (x)| 1 f(x) = 6-2 F + y, 


2) 13 
1 31 | 6) 
= . => — 2 —————_ 4 — aialiel 
F g(x)" y + f(z), yy 240 © I+ eo 4ag0° & Smee 
The derivative y’ can be found from the relation 
i a 19 
— ie es ily He ste one. Cae bes 
hy’ =ypo'F 1g hF g HOY t+ arg Ko'y ab (7) 


which is easily derived from (4) by using the central difference formula for F’ 
and eq. (6) for the difference correction y. According to (6) and (7) the initial 
values of the first and second sums of F are then given by® 


(8) 

61F, = 17 + het You ae as oe ea es 

cD ys Ce he ee ee ne 
Let us now first consider the case when the interval h has been choosen so 
small that the difference correction y can be entirely neglected in the whole 


range of integration. The numerical integration of (5) can then be carried out 
by using the recurrence rule indicated by the arrow in the scheme below: 


-2 

Xo Yo 0" Fo é-1 Fy Fo 

Ly YW 6? Fy 

Lyi Mey hoe OF” Be Fy (9) 
a Fri 


Soci hn eae) Tee Migs 1 


| 


If the value of y», has been recorded, it is possible to compute F,, from eq. 
(6”) and to find 6-!F,4, and 6-?Fy+41 by successive summations: 


Oba Oa Py, CO ln =O ln FO Enh. (10) 


7 W. E. Mitne, Amer. Math. Monthly 49, 96 (1942); L. Fox and E. T. Goopwiy, Proc. 
Camb. Phil. Soc. 45, 373 (1949). 

8 Compare B. Noumeroy, Publ. Observatoire Astrophys. Cent. Russie, vol. ii, Moscou 
1923; Mon. Not. Roy. Astr. Soc. 84, 592 (1924). 

9 In the following we will use the abbreviations y, = y (x) + rh) ete. 
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Then yn+1 can be obtained according to (6’) by an addition and a division, and 
the process can be repeated for the next step. If the functions f(x) and g (a) 
were tabulated in advance, we have found it possible to carry through this 
process by the aid of only ordinary desk machines with nine figures in the — 
calculation with a rapidity of 60—80 steps an hour. 

In comparison to the methods proposed by Hartres*t, by Hartrer, Kronic 
and PETERSEN”, or by BuckrncHam®, we note that here no Taylor expansion 
is used in the start and that further no estimations are necessary. 

In many problems, in which the integration has to be extended over a wide 
range demanding an enormous number of steps, this manner of proceeding can 
be rather laborious. In such cases it seems better to use a larger interval h 
and to take the difference correction y into full account. In considering y it 
is suitable to divide the integration into three parts: the start, the marching 
process, and the aftercorrection. 

In the start region we will determine the correction y by using a method of 
successive approximations. We let an upper index indicate the order of ap- 
proximation. In a zero-order approximation we put y°=0, and some values 
of y® and F® around the starting point are then computed according to (6) 
and (9). From the values of F® and its differences a first-order approximation 
y® of the difference correction can be determined, and the start is then re- 
peated until two consecutive series of y” agree. In most cases the first ap- 
proximation appears to have a sufficient degree of accuracy. 

When a sufficient number of values of F has been found, we will proceed 
in another way in order to spare work. Let us assume that the values of y 
and F up to the point «= 2% +nh have been computed, and that the cor- 
responding difference scheme of F has been recorded. If V denotes the back- 
ward difference operator V = 1— H™!, the quantities in the last backward line 
m in this scheme may be denoted by Fy, V Fn, V?2Fn, V32Fn,... They may 
be used for extrapolating the difference correction, for from the symbolic rela- 
tio #=(1—V)!=1+V+V2+--- and (2) we easily derive the extrapola- 
tion formulas 


02 Fait = V2 Fn + 2V2 Fn + 8V4 Fn + 4VERn + -*, 


OOF, = VAP, + BVEP, 6 Vokes ELV IR, 2 5, co 
and then ynx+i is formed according to (6’’’). We note that the elementary series 
(11) are suitable for calculations by means of ordinary desk machines.!° By 
using (6), (10), and (11), it is now possible to carry through the integration in 
a quite straightforward way, and we will therefore characterize this part of the 
calculation as a “marching process”, which gives a preliminary solution Ypret. 

Afterwards, when the marching process has been completed for the whole 
range of integration (including some additional points), it is possible to improve 


10 In order to be fit for use the series (11) must be rapidly convergent, which can be 
obtained by choosing the interval h not too large. However, in general it is convenient to 
calculate the differences of # only up to the order p where the irregularities, due to the 
effect of rounding-off errors, become to be of about the same magnitude as the differences 
themselves; then the formulas (11) correspond to an extrapolation of F and its differences 
by a polynomial of the degree p. 
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the accuracy of the solution by utilizing the actual values of 6? F,,, 64Fn, 
found on later stages of the calculation. The final solution is then determined 
by the formula 


Yn, tinal = Yn, prel + h? (Yn, actual — Yn, extrapolated). (12) 


The idea of the method is that the “‘aftercorrection” in (12) shall be so small 
that, in reality, it will be only a check of the extrapolations (11), possibly in- 
fluencing the last figure in y. The interval h can therefore not be choosen 
too large. 

If the difference correction is taken into full account in this way, we have 
found it possible to carry out the numerical integration by means of only 
ordinary desk machines with nine figures in the calculation with a rapidity of 
20—30 steps an hour. 

Finally, in comparison to the methods proposed by MiLnE? and by Fox and 
GoopwINn’, we note that the solution y is here obtained directly in summed 
form, which is an advantage; see the discussion by Jackson. 


Ill. Calculation of eigenfunctions and asymptotic phases 


The radial wave equation of the scattering problem connected with the cen- 
tral field V(r) can be written in the form 


Ch re bia i 


ie eae rma | 


wy! ss 
v(pP = 0, (13) 


where v(r) represents the potential and k& is the wave number: 


822m 8ma?m 


Liheke aa © 


j2 = V(r). (14) 


Following HuLrHEN we will assume that v(r) is at least sectionally continuous 
and bounded in the whole domain except at r = 0, where it may be singular 
as c,/2” —c,/x2, p< 2, c, arbitrary, and cg 20. For the sake of simplicity we 
will further assume that rv(r) approaches zero as r tends to infinity, but we 
note that the following discussion without difficulties can be generalized to in- 
clude even the Coulomb case. 

The solution with the initial condition ¥Y (0) = 0, Y’(0) = 1, has the asymp- 
totic form 


W (r) ~ A sin ( ) (15) 


where 7 is the asymptotic phase giving the partial cross-section connected with 
the /-scattering. 

Starting from the initial condition, we can integrate the wave equation (13) 
numerically be means of the method described in Sec. II. The integration has 
to be performed up to the point r= 7,, where v(r) begins to be negligible in 
comparison to k?, and for r=vr, the solution has then the form 
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YW (r) = const. {sin 71 Gs (7) + cos m Gy (r)} (16) 
with ° . — 
G, = (737) Jisy (kr), Gy =(—1) (7) J-1-4 (kn), (17) 


Here’ J» (x) are the Bessel functions of the first kind. The asymptotic phase 
ym can then be determined from the relation 


PG ai re) 


= 1 


where the derivative Y’ may be computed from (4) or from the linear alge- 


lod 


braic equation (7). A good check on the calculations is obtained by computing 


the constant 7 for different values of r. 

In order to investigate the propagation of errors in our integration scheme (9), 
we have also considered another way of integration. Following Mune" we 
put Y =o sin y with the auxiliary condition 7’ = «w~*, « = arbitrary constant, 
which for @ gives the equation 


d? w L(E+ 1) 
ye was 
dr? a tf r 


ae v(nfo— =.=, (19) 


which again can be numerically integrated by the Gaussian formula (3). This 
equation has previously been used for investigating the asymptotic behaviour 
by Rogansky and WeErTzEL?2, by Morsse?, by Hartree, Kronic and PETERSEN?, 
and by Arnot. We have integrated eq. (19) enwards, starting from the initial 
condition w (co) = («/k)?, w’(co) = 0. We have then connected the solution 
with the wave function ¥Y, found by starting from r = 0, in an arbitrary point 
r=re by the formula 


; cot x. (20) 


The asymptotic phase 7 can then be found from the relation 


] 2 , y 7 a . 
Tier = arccot fe (5 =\l kre + xf (o2—*)ar (21) 


1g 


c 


The formulas (20) and (21) were first given by Harrrer, Kronic, and PrrersEn.? 
We note that the result must be independent of 7- and of the value of «; for 
the sake of simplicity we have choosen « = k. We have compared the different 
values of 2 obtained according to (18) or to (21) in some typical cases, see 
Table I, and, since the agreement was excellent, we have concluded that the 
propagation of error in the scheme (9) could be entirely neglected. For values 


1 W. E. Mitne, Phys. Rev. 35, 863 (1930); H. A. Witson, Phys. Rev. 35, 948 (1930); 
E. L. Hitt, Phys. Rev. 38, 1258 (1931); L. A. Youne, Phys. Rev. 38, 1612 (1981). 
E. Mavetune, Zeits. f. Physik 67, 516 (1931) has treated the special case « = 1. 

2 F, L. Arnot, Proc, Camb. Phil. Soc. 32, 161 (1936). 
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of k, which are not too large, we have therefore found it convenient to use 
formula (18) alone. It is, of course, much simpler to integrate the linear 
equation (13) by means of (6) than to integrate the non-linear equation (19) 
by using (3). 


IV. Applications to the Yukawa potential (S-phases). Comparison with the 
results obtained according to Hulthén’s variational principle 


According to HuttuEn® and others! 1516 the asymptotic phases or the 
solution of the scattering problem can be determined from a variational prin- 
ciple by using a series of approximate wave functions Y,, Y,, WY; ... with an 
increasing number of parameters c,. HunTHén considered the integral 


L= [PLP dr, (22) 
0 
for all functions Y having the boundary conditions 
: la 
W (0) = 0, ¥ (r)~ sin (kr — + n) (23) 


with different phases 7. Varying the function Y around the exact solution VY, 
HULTHEN obtained the formula 


dL +kbn=0. (24) 


If Y is varied in such a way that the additional condition ZL = 0 is always 
fulfilled (6Z = 0), this gives 67 = 0, which means that the asymptotic phase 
has a stationary character. Since, in general, the “eigenphase” is not an ex- 
tremum, the variational principle must in this case be treated with special 
caution, and HuitHen has worked out two criterions for the accuracy of the 
results obtained. 

HULTHEN’s original formulation of the variational principle was rather com- 
plicated in the applications, since the equations for tg or coty appeared to 
be of the second degree. Recently he has therefore developed a new modifica- 
tion!” of his method for calculating S-phases, which is considerably simpler. The 
new method consists of two steps, of which the first is again based on a varia- 
tional principle. By dropping the condition L = 0, HuLTHEN obtains now linear 


13 LL. Huirunn, K. Fysiogr. Sallsk. Lund Férhandl. 14, no. 21 (1944); ibid. 14, no. 8 (1944); 
Cc. R. X™° Congrés des Math. Scan. 1946, Copenhagen 1947, p. 201; compare also Arkiv f. 
mat., fysik o. astr. 29 B, no. 1 (1942). 

14. A. Hytieraas, Fra Fysikkens Verden (Nie~s Bour issue), Oslo 1945; A. Pats, Proc. 
Camb. Phil. Soc. 42, 45 (1946); N. Zenon, K. Fysiogr. Sallsk. Lund Forhandl., N.F. 58, 
no. 10 (1947); I. E. Tamm, J. Exp. Theor. Phys. 18, 337 (1948); ibid. 19, 74 (1949). 

15 W. Koun, Phys. Rev. 74, 1763 (1948); 8S. S. Huane, Phys. Rev. 76, 477 (1949); ibid. 
76, 867 (1949). 

16 J, ScHwINGER, Phys. Rev. 72, 742 (1947); “Lectures on Nuclear Physics’, Harvard 
University 1947; J. M. Brarr—J. D. Jackson, Phys. Rev. 76, 18 (1949); B. A. Lippmann— 
J. Scuwincer, ibid. 79, 469 (1950); B. A. Lippmann, ibid. 79, 481 (1950). 

17 .. Hunrutn, Arkiv f. mat., fysik o. astr. 35 a, no. 25 (1948). 
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equations for the parameters ¢,. The first step gives only preliminary values 
y’ and L’ for the phase and the integral (22), respectively, but the phase 
value can be essentially improved in the second step. ,Puttmg 079 = 7’ —7'' 
and 6L = L'—L" = L’ into formula (24), we get immediately the relation 


1 
ie = 7 + ZU (25) 


iad 


where 7’ is the improved phase value.!8 

In this modified form of HuLrHEn’s variational principle with its two steps, 
it is not so easy to get an analytical test of the accuracy of the final result 
for the asymptotic phase, even if HULTHEN’s two criterions give a rough estima- 
tion of the accuracy of the wave functions. In such a case, it may therefore 
be of some interest to study the final accuracy in some typical results by 
means of numerical methods. 

In order to carry out a comparison with HuLTHEN’s results’, we have treated 
eq. (13) with the Yukawa potential v(r) = 1.5e7"/r, belonging to the deuteron 
problem, for 1/=0 and for k= 0.1, 0.2,...,1.0. The numerical integration 
was carried out according to Sec. II with eight or nine figures in the calcula- 
tion and the intervals h = 0.1 and h = 0.2. By using ordinary desk machines 
(Facit ESA), it was possible for us to record a wave function for a given k- 
value and to find the corresponding phases within 8—9 hours, including the 
start and the check of the integration. The results are contained in Tables 
I—III and in Figs. 1 and 2. 


Table I 


A study of the numerical accuracy and of the propagation of errors in the 
integration scheme (9), obtained by using outward and inward solutions con- 
nected in the point r= 7; see formula (21) 


. 


v(r)=1.5e7"/r, 1= 0, b= 0.2, I= kf (w-2—1)dr 


ie cot (X¥,—kr,) I cot 
6 0.671 1752 0.000 1511 0.670 9561 
8 0.670 9801 165 0.670 9562 
10 0.670 9590 19 0.670 9563 
12 0.670 9566 2 0.670 9563 
| 14 0.670 9563 0 0.670 9563 


Table II contains a comparison between Hunrnin’s phases!? and the values 
computed by us. It is perhaps somewhat surprising that so extremely good 
‘phase values can be obtained by using Hutrutn’s modified principle when 


** It may be noted that, by putting 7’ = 0, Y= G, into (25), one gets the well known 
formula for 7 in Born’s approximation. 
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Fig. 1. 


only one, two, or three variational parameters are used. We note that the 
first-order improvement (25) is essential for this result. 

Finally we have investigated how HvuLruin’s approximate wave functions 
tend against the accurate function obtained by numerical integration. Hut- 
THEN’s functions are of the form 


Yn(r) = sin (kr + 4) — Yn (r)- sin 7, 
(26) 


n 
dn = ase (ite) op 0m" 
y=1 


where yn(r) gives the deviation from the ‘‘asymptotic” form in the neigh- 
bourhood of r=0. Figures 1 and 2 contain for k = 0.1 and k = 1.0, respec- 
tively, the function y(r) found from the numerically calculated wave function 
V(r) by a normalization to the amplitude A= 1 for large values of r. In 
another scale the same figures contain also the differences yp(r)—y(r) for 
n = 1, 2, 3, ie. the errors in HuLTHEN’s approximations. It appears that Hut- 
THEN’s functions y,(r) give good approximations of y(r), and that the con- 
vergency is almost uniform in the whole range for r. A good wave function 
W(r) can therefore be obtained by using (26), the constants c, found from 
the variational principle in the first step, and the phase value 7 found from 
(25) in the second step. 
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Table II 


Comparison between Hutrutn’s phase values!? and the values obtained by 
using numerical integration according to Sec. II and formula (18). = number 
of variational parameters c, 


Hulthén’s values 
k n cot 7) cot 1) cot 7) difference 
variational first-order numerical (2) — (3) 
principle improvement integration 
(1) (2) (3) (4) 
0.1 1 1.021 92 0.984 50 0.000 90 
2 0.999 36 0.983 67 27 
3 0.990 84 0.983 46 6 
== 0.983 395 \ 
0.2 1 0.698 03 0.671 87 0.000 91 
2 0.682 49 0.671 19 23 
3 0.676 55 0.671 02 6 
— 0.670 956 
0.3 1 0.665 21 0.639 57 0.001 O1 
2 0.650 14 0.638 81 25 
Be 0.644 45 0.638 63 7 
= 0.638 563 
0.4 ] 0.698 73 0.670 87 0.001 12 
2 0.682 42 0.669 99 24 
si 0.676 30 0.669 82 7 
— 0.669 753 
0.5 1 0.754 17 0.722 74 0.001 15 H 
2 0.735 96 0.721 80 21 | 
ap 0.729 15 0.721 64 5 
— 0.721 594 
0.6 1 0.817 52 0.781 52 0.001 06 
2 0.796 89 ‘0.780 60 14 
3 0.789 19 0.780 47 1 
== 0.780 455 
0.7 1 0.883 73 0.842 04 0.000 87 
By 0.860 21 0.841 25 8 
3 0.850 89 0.841 17 0 
ee 0.841 171 
0.8 1 0.950 99 0.902 21 0.000 58 
2 0.924 19 0.901 65 2 
3 0.91418 0.901 63 0 
= 0.901 629 
0.9 1 1.019 14 0.961 20 0.000 24 
2 0.990 12 0.960 96 0 
3* 0.973 53 0.960 99 3 
= 0.960 961 
1.0 1 1.089 28 1.018 71 —0.000 14 
2 0.992 61 1.018 97 =+-12 
3 1.033 87 1.018 92 7 
= 1.018 854 


* Private communication from Prof. L. Hulthén. 
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Iie, Pe 


A comparison between Table III, which contains HULTHEN’s criterion ratios 
0; and gg, and Figs. 1 and 2 shows that the quantities |0,—1]| and |e,—1| 
can roughly be taken as measures of the errors yn — y. 


Table III 


HULTHEN’s criterion ratios 0; and 0, for k= 0.1 and k= 1.0. These quantities 
are equal to 1 for the exact solution. m = the number of variational parameters 


| E— Ost B10 

| n | Q,—1 Q,—1 | Q,—1 | Q,—1 

| i 

| | | 

1 | —0.00 55 —0.00 89 | —0.04 26 | —0.00 64 

| 2 —0.00 23 | +0.00 48 -+-0.01 59 | +0.01 54 
3 —0.00 11 | —0.00 29 —0.00 90 | +0.00 34 


V. Comments 


We have found that, by means of the methods developed in Secs. II and 
III, it is possible to calculate a wave function and the corresponding asymp- 
totic phase with an accuracy of 3—4 figures within a few hours and with an 
accuracy of 6—7 figures within 8—9 hours by using only ordinary desk ma- 
chines. 
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Our study has further shown that Hutruién’s modified form of the varia- 
tional principle!’ with the improvement formula (25) can give excellent phase — 
values even when only a few variational parameters are used. For the Yukawa. 
potential and / = 0, HutrH&n has given approximate wave functions which are — 
well fitted for the purpose; for other potentials the construction of these trial 
functions seems to be one of the essential problems. 


We are greatly indebted to Professor Lamek Hulthén, Stockholm, for placing 
some of his unpublished numerical results at our disposal, and for some valu- 
able comments on Part IV. 


Institute of Mechanics and Mathematical Physics, University of Uppsala, 
Uppsala, Sweden. 


Tryckt den 23 april 1951 


Uppsala 1951. Almqvist & Wiksells Boktryckeri AB 
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